#include "Nurbs.h"

/********************************************************************************************************
Constructor & Destructor
********************************************************************************************************/

//Default Constructor
Nurbs::Nurbs(){

}

Nurbs::Nurbs(float* vertex, int nbPoints, int degree, float* weight, float step){

	m_lineClose = isLineClose(vertex, nbPoints);
	m_nbPoints = 0;
	computePassagePoints(vertex, nbPoints, degree, weight, step);

}


Nurbs::~Nurbs(){

}


/********************************************************************************************************
Compute Nurbs Points
********************************************************************************************************/

void Nurbs::computePassagePoints(float* vertex, int nbPoints, int degree, float* weight, float step){


	int* knots = computeIntervals(nbPoints, degree);
	int nbKnots = nbPoints + degree;
	
	std::vector<float> newPoints = std::vector<float>();
	float* point;

	//First Point
	newPoints.push_back(vertex[0]);
	newPoints.push_back(vertex[1]);
	newPoints.push_back(vertex[2]);
	m_nbPoints++;

	for(float t = (float)knots[0]; t < (float)knots[nbKnots-1]; t += step)
    {
     
		if(t != 0 )
		{
		  point = evaluate(vertex, nbPoints, knots, degree, weight, t);
		  newPoints.push_back(point[0]);
		  newPoints.push_back(point[1]);
		  newPoints.push_back(point[2]);
		  m_nbPoints++;
		}
   
    }
	
	//Last Point
	newPoints.push_back(vertex[nbPoints * 3 - 3]);
	newPoints.push_back(vertex[nbPoints * 3 - 2]);
	newPoints.push_back(vertex[nbPoints * 3 - 1]);
	m_nbPoints++;


	m_vertex = new float[m_nbPoints * 3];
	for(int i = 0; i < m_nbPoints; i++){
		m_vertex[(i * 3)] = newPoints[(i * 3)];
		m_vertex[(i * 3) + 1] = newPoints[(i * 3) + 1];
		m_vertex[(i * 3) + 2] = newPoints[(i * 3) + 2];
	}

	
}


//Compute Knots
int* Nurbs::computeIntervals(int nbPoints, int degree)
{	
  int *intervale = new int[nbPoints + degree];
  
  for (int i = 0; i < nbPoints + degree; i++)
  {
    if (i < degree)
      intervale[i] = 0;
    else
    {
      if (degree <= i && i <= nbPoints)
        intervale[i] = i - degree + 1;
      else if (i > nbPoints)
        intervale[i] = nbPoints - degree + 2;
    }
  }
  
  return intervale;  
}


//Evaluate point
float* Nurbs::evaluate(float* vertex, int nbPoints, int* knots, int degree, float* weight, float t){

	float* point = new float[3];
	point[0] = 0;
	point[1] = 0;
	point[2] = 0;
	
	float denum = 0;
	int i;
	for(i = 0; i < nbPoints; i++){
	  denum += computeValue(knots, t, i, degree+1) * weight[i];
	}
  
	for(i = 0; i < nbPoints; i++){
		float compute = computeValue(knots, t, i, degree+1);
	  point[0] = point[0] + vertex[(i*3)] * compute * weight[i];
	  point[1] = point[1] + vertex[(i*3) + 1] * compute * weight[i];
	  point[2] = point[2] + vertex[(i*3) + 2] * compute * weight[i];
	}
  
	if(denum != 0.0f){
		point[0] /= denum;
		point[1] /= denum;
		point[2] /= denum;
	}

	return point;
}

/*
float Nurbs::computeValue(int* knots, float u, int i, int p)
{
	float basis_i;
	float basis_i_1;
	float temp1, temp2;
	
	//Stop recursivity
	if(p == 0)
    {
      //Check if near a control Point
      if(knots[i] <= u && u <= knots[i+1])
		return 1;
      else 
		return 0;
    }

	else
    {
      
		if((knots[i + p] == knots[i]) && (knots[ i + p +1] == knots[i + 1]))
			return 0;
		else
		{
		    if((knots[i + p] == knots[i]) || (basis_i = computeValue(knots, u, i, p-1)) == 0)
				temp1 = 0;
			else
				temp1 = ((u - knots[i]) / (knots[i + p] - knots[i])) * basis_i;
			if((knots[ i + p +1] == knots[i + 1] || (basis_i_1 = computeValue(knots, u, i + 1, p - 1)) == 0))
				temp2 = 0;
			else
				temp2 = ((knots[i + p + 1] - u) / (knots[ i + p +1] - knots[i + 1])) * basis_i_1;
		}
    }

  return (temp1 + temp2);

}
*/

float Nurbs::computeValue(int* intervals, float value, int step, int degree){
    float result = 0;
    
    if (degree == 1)
    {
      
      if (intervals[step] <= value && value < intervals[step + 1])
        result =  1;
      else
        result = 0;
    }

    else
    {
      if (intervals[step + degree - 1] == intervals[step] && intervals[step + degree] == intervals[step + 1])
        result = 0;
      else if (intervals[step + degree - 1] == intervals[step])
        result = (intervals[step + degree] - value) / (intervals[step + degree] - intervals[step + 1]) * computeValue(intervals, value, step + 1, degree - 1);
      else if (intervals[step + degree] == intervals[step + 1])
        result = (value - intervals[step]) / (intervals[step + degree - 1] - intervals[step]) * computeValue(intervals, value, step, degree - 1);
      else
      {
        result = (value - intervals[step]) / (intervals[step + degree - 1] - intervals[step]) * computeValue(intervals, value, step, degree - 1) + (intervals[step + degree] - value) / (intervals[step + degree] - intervals[step + 1]) * computeValue(intervals, value, step + 1, degree - 1);
      }
    
    }


    return result;
}
